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Abstract 



The theory of quasi-arithmetic means is a powerful tool in the study of 
covariance functions across space-time. In the present study we use quasi- 
arithmetic functionals to make inferences about the permissibility of aver- 
ages of functions that arc not, in general, permissible covariance functions. 
This is the case, e.g., of the geometric and harmonic averages, for which we 
obtain permissibility criteria. Also, some important inequalities involving 
covariance functions and preference relations as well as algebraic properties 
can be derived by means of the proposed approach. In particular, we show 
that quasi-arithmetic covarianccs allow for ordering and preference relations, 
for a Jensen-type inequality and for a minimal and maximal element of their 
class. The general results shown in this paper are then applied to study of 
spatial and spatiotemporal random fields. In particular, we discuss the rep- 
resentation and smoothness properties of a weakly stationary random field 
with a quasi-arithmetic covariance function. Also, we show that the gener- 
ator of the quasi-arithmetic means can be used as a link function in order 
to build a space-time nonseparable structure starting from the spatial and 
temporal margins, a procedure that is technically sound for those working 
with copulas. Several examples of new families of stationary covariances 
obtainable with this procedure are shown. Finally, we use quasi-arithmetic 
functionals to generalise existing results concerning the construction of non- 
stationary spatial covariances and discuss the applicability and limits of this 
generalisation. 

KeyWords: Completely monotone functions, Jensen's type inequalities, 
Nonseparability, Nonstationarity, Quasi-arithmetic functionals. Random Fields, 
Smoothness, Space-time covariances. 
Subj-class: Probability; Statistics. 
MSG class: 60G60; 86A32. 
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1 Introduction 



The importance of quasi-arithmetic means has been well understood at least 
since the 1930s, and a number of writers have since contributed to their 
characterisation and to the study of their properties. In particular, Kol- 
mogorov (1930) and Nagumo (1930) derived, independently of each other, 
necessary and sufficient conditions for the quasi-arithmeticity of the mean, 
that is, for the existence of a continuous strictly monotonic function / 
such that, for xi, . . . , ,x„ in some real interval, the function (xi, . . . , Xn) ^ 
Mn{xi, ...,Xn) = f^^Ti YJifi^i)) is & mean. 

Using this result, they partially modified the classical Cauchy (1821) 
internality and Chisini's (1929) invariance properties. As pointed out by 
Marichal (2000), the Kolmogorov reflexive property is equivalent to the 
Cauchy internality, and both are accepted by statisticians as requisites for 
means. 

Early works on the concept of mean include Bonferroni (1926), de Finetti 
(1931), Gini (1938), Dodd (1940) and Aczel (1948). More recent contribu- 
tions are the works of Wimp (1986), Hutnik (2006), Matkowski (1999;2002) 
Jarczyk and Matkowski (2000), J. Marichal (2000), Daroczy and Hajdu 
(2005), and Abrahamovic et al. (2006). Quasi-arithmetic means, in par- 
ticular, have been applied in several disciplines. Their functional form has 
been used in the theory of copulas under the name of Archimedean copulas 
(Genest and MacKay, 1986) and a rich literature can be found under this 
name. In the theory of aggregation operators and fuzzy measures, a growing 
literature related to the use of quasi-arithmetics includes the works of Prank 
(1979), Hajek (1998), Kolcsarova (2001), Klement et al. (1999), Grabish 
(1995) and Calvo and Mcsiar (2001). 

Despite the extensive quasi-arithmetic means literature, to the best of 
our knowledge, there is no published work relating quasi-arithmetic means 
with covariance functions, whose properties have been extensively studied 
both in mathematical analysis and statistical fields. In particular, the study 
of covariance functions is intimately related to that of positive definite func- 
tions, the latter being the subject of a considerable literature in a variety 
of fields, such as mathematical analysis, abelian semigroup theory, spatial 
statistics and geostatistics. For basic facts about positive definite functions, 
we refer the interested reader to Berg and Forst (1975) and to Berg et al. 
(1984). The importance of positive definite functions in the determination of 
permissible covariance functions (ordinary and generalised) in spatial statis- 
tics was studied in detail by Christakos (1984). Subsequent considerations, 
in a spatial and a spatiotemporal context, include the works of Christakos 
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(1990, 1991, 1992), Sasvari (1994) and Gneiting (1997), among others. 

Fundamental properties of covariance functions may be inferred by study- 
ing collections of them considered as convex cones closed in the topology of 
point-wise convergence. In this paper, quasi-arithmetic averages and pos- 
itive definite functions are combined to gain valuable insight concerning 
certain covariance properties. Also, we apply the general results obtained 
by our analysis based on the concept of quasi-arithmcticity to build new 
classes of stationary and nonstationary space-time covariance functions. In 
such a context, we seek answers to the following questions: 

(1) Consider an arbitrary number n G N of covariance functions, not neces- 

sarily defined in the same space. Their arithmetic average and product 
{i.e., the geometric product raised to the n-th power) are valid covari- 
ance functions, and so is the fc-th power average of covariance functions 
{k is a natural positive number). But, what about other types of av- 
erages? Since quasi-arithmetic means constitute a general group that 
includes the arithmetic, geometric, power and logarithmic means as 
special cases, it seems natural to use quasi-arithmetic representations 
in order to derive positive definiteness conditions for quasi-arithmetic 
averages of covariance functions. 

(2) Can we use the properties of quasi-arithmetic means to establish impor- 

tant inequalities, ordering and preference relations and minimal and 
maximal elements within the class of covariance functions? 

(3) Is it possible to find a class of link functions that, when applied to the 

k covariances, can generate valid nonseparable covariances? If this 
is the case, other potentially desirable properties should be examined, 
e.g., this approach should be as general as possible and include famous 
constructions and separability as special cases; and it should preserve 
the margins. 

In view of the above considerations, the paper is organised as follows: 
Section 2 discusses the background, notation and proposed methodology; it 
also provides a very brief introduction to positive definite functions. In Sec- 
tion 3, the main theoretical results are presented. In particular, we propose 
permissibility criteria for the quasi-arithmetic mean of covariance functions. 
We derive important covariance inequalities of Jensen's type as well as order- 
ing and preference relations between covariance functions. Finally, minimal 
and maximal elements of the quasi-arithmetic covariance class are identi- 
fied, and an associativity property of this class is provided. In Section 4, 
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we apply our results to the construction of new families of space-time non- 
separable stationary covariance functions. Also, we extend Stein (2005b) 
result to a more general class of spatial nonstationary covariance functions. 
Several examples of stationary and nonstationary covariances are proposed 
and their mathematical properties discussed. Finally, we study the proper- 
ties of quasi-arithmetic random fields (defined later in the paper) in terms 
of mean square differentiability and variance. Section 5 concludes with a 
critical discussion of the preceding analysis. All the proofs of the original 
results derived in this paper can be found in the Appendix. 

2 Background and Methodology 

2.1 Covariance functions: characterisation and basic prop- 
erties 

In this section, we shall assume the real mapping C:XxXCM"'x]R^^]R 
to be continuous and Lebesgue measurable on the domain X x X, where X 
can be either a compact space or the entire d-dimensional Euclidean space, 
d € N. C is a covariance function if and only if it is positive definite, that is 

n n 

5^^CiC,C(xi,x,) >0 (1) 

i=l i=l 

for any finite set of real coefficients {cj}"^]^, and for xi, . . . ,x„ € X. Chris- 
takos (1984) calls the covariance condition ([T]) permissibility, and throughout 
the paper we shall use both this term and that of positive definiteness to 
characterise valid covariance functions. A subclass of positive definite func- 
tions, called stationary, is obtained if 

C(xj,Xj) := Co(x), 

with X = Xj — Xj and Co : X ^ M positive definite and such that Co(0) < oo. 

As shown by Bochner (1933), condition ([1]) is then equivalent to the re- 
quirement that Co is the Fourier transform ^ of a positive bounded measure 
Co with support in W^, that is 

Co(x):=.F[Co](x)= / e^'^'''dCo(c^). 

Additionally, if Co is absolutely continuous with respect to the Lebesgue 
measure (ensured if Cq G Li(X)), then the expression above can be written 
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as a function of dCo(w) = co(a;)da;, where cq is called the spectral density 
of Co- For a detailed mathematical discussion of the Fourier representation 
in a spatial-temporal statistics context, we refer the interested reader to 
volumes by Yaglom (1987) and Christakos (1992). 

In what follows, we shall consider some interesting restrictions on the 
general class of stationary covariance functions. But first, we need to intro- 
duce some standard notation for arbitrary partitions and operations between 
vectors. In order to manipulate arbitrary decompositions of nonnegative in- 
teger numbers, let us consider the set exp(Z+) = U Z_|_ U U U . . . 
(disjoint union). An element d of exp(Z+) can be expressed either as d = 
or as d = (di,d2, ■ ■ ■ ,dn) if d G Z" with n > 1. In the latter case we denote 
by n(d) = n the dimension of d and |d| = Y17=i length of d. Both 

values are taken to be whenever d = 0. For d,d' G exp(Z_i_) we say that 
d < d' if and only if n(d) = n(d') and di < d[ for alH = 1, 2, . . . ,n(d). 
Usual vector operations are possible only between elements of the same di- 
mension. Vectors with all components equal are denoted in bold symbols, 
such as 0, 1. 

Now, the restrictions considered in this paper are of the following type: 

1. The isotropic case: 

Co(x) := Ci(||x||), xgX, (2) 

that is, C is said to be represented by the positive definite function Ci : 
R_(- M, that is rotation-translation invariant (or radially symmetric) , 
and where ||.|| denotes the Euclidean norm. This is the most popular 
case in spatial and spatiotcmporal statistics (Matern, 1960; Yaglom, 
1987, Christakos and Papanicolau, 2000). 

2. The component-wise isotropic case: Let us consider the ti-dimensional 
space W^, and let d be an element of exp(Z_|_) such that |d| = d and 

1 < d. Thus, one can create opportune partitions of the spatial lag 
vector X G M*^ in the following way. If d = {di,d2, ■ ■ ■ ,dn) and x G R'' 
we can always write 

X = (xi, X2, . . . , x„) G M*^! X R''^ X ■ ■ ■ X W^" 

so that: 

(i)- Co(x) = Co(k) for any x, k G M*^ if and only if ||xj|| = ||kj|| for all 
i = 1, 2, . . . , n. 
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(ii)- The resulting covariance admits the representation 
Co(x) := Ci(||xi||,...,||x„||) 

/"OO /'OO 

= /.../ J[^dA\\^i\\n)dF{n,....rn) (3) 
Jo Jo 

with nd{t) = T{d/2) {j) J(d-2)/2{t), Jd{-) denoting the Bessel func- 
tion of the first kind of order d (Abrahamowitz and Stegun, 1965), 
F a n-variate distribution function and Ci : M" — > M positive defi- 
nite. Thus, equation is a special case of ^ (jd| = d and d := d, 
a scalar) and its corresponding integral representation can be read- 
ily obtained. The special case d = 1 and |d| = d is particularly 
interesting in the subsequent sections of this paper, as the function 
{xi, . . . ,Xn) I— >■ Ci(|xi|, . . . , |xn|), with Ci positive definite, does not 
depend on the Euclidean norm, but on the Manhattan or city block 
distance, with important implications in spatial and spatiotemporal 
statistics as pointed out by Christakos (2000) and Banerjee (2004). 

It is worth noticing that covariance functions of the type ([3]) have a property 
called reflection symmetry (Lu and Zimmerman, 2005), or full symmetry 
(Christakos and Hristopulos, 1998). This means that C(xi, . . . , Xj, . . . , x„) = 
C(xi, . . . , -Xj, . . . ,x„) = ... = C(-xi, . . . , -Xj, . . . , -x„). Whenever no 
confusion arises, in the remainder of the paper we shall drop the under- and 
super-script denoting a stationary or stationary and isotropic covariance 
function, respectively. 

As far as the basic properties of covariance functions are concerned, 
assume that Cj : M'^* — > M+ (i = 1, . . . , n) are positive, continuous and inte- 
grable stationary covariance functions with di G Z+ and let d = (di, . . . ,dn) 
such that |d| = d. It is well known that some mean operators preserve the 
permissibility of the resulting structure. In particular, if we assume, without 
loss of generality, that the 9i are nonnegative weights summing up to one, 
the following are then permissible on W^: 

1 The arithmetic average 

n 

Ca(x) = ^e,C,(x,), 

i=l 

2 the non weighted geometric average up to a power n, 

n 

Cg(x) = J]Q(x,), 
1=1 
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3 the A;-power average {k G '^+), up to a power A;, 

n 
i=l 

4 Scale and power mixtures of covariance functions (Christakos, 1992), 

C(x)= / C(x;e)d//(e), 

for jj. a positive measure and G 9 C W, p G N. 

2.2 The methodology. Quasi-arithmetic multivariate com- 
positions 

Quasi-arithmetic averages have been extensively treated by Hardy et al. 
(1934). Our methodology generalises the concept of quasi-arithmetic aver- 
ages and introduces new formalisms and notation, since our aim is to find a 
class of compositions of covariance functions that satisfies desirable proper- 
ties. 

Let ^> be the class of real-valued functions (p defined in some domain 
D{ip) C M, admitting a proper inverse (p~^ defined in D{ip~^) C M and 
such that ip{ip~^{t)) = t for all t G D{ip~^). Also, let $c and $cm be 
the subclasses of <1> obtained by restricting ip to be, respectively, convex or 
completely monotone on the positive real line. Let us call a quasi- arithmetic 
class of Junctionals the class 

£1 := 1^ : D{ip-'^) x ■ ■ • x D{ip-^) ^ M : ^(u) = if eiip-^{ui)^ , <^ G $ 

. 

where 9i are nonnegativc weights and u = (ui, . . . ,Un)', for n > 2 positive 
integer. Also, we shall call Oc and £2cm the subclasses of when restricting 
ip to belong, respectively, to $c and $cm- 

If V' G 0, then we should write ^p^ as the function such that: for any 

nonnegative vector u, '^^{vl) = (p^ (j2i=i ^iV^^ {ui)^ . For ease of notation, 
we simply write ip instead of ip^, whenever no confusion arises. 

Next, we introduce a new class of functionals that will be used extensively 
throughout the paper. 

Definition 1. Quasi- arithmetic compositions // fi : M+ such 

that U"/i(M''*) C D{ip~^) for some G the quasi-arithmetic composition 
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Table 1: Examples of quasi-arithmetic compositions for some possible 
choices of the generating function (/? G <I>. 



fit) 




Q^(/l,/2)H 


Remarks 


exp(— t) 


-Int 




fi : M'*' ^ [0, oo) 

J I L ■ / 

InO = —oo 
exp(— oo) = 

Efi = i 


1/t 


1/t 




fi-.R""^ ^ [0,cx)) 
1/0 = oo, l/oo = 

0/0 = 


M(l -t/M)+ 


M{l-t/M)+ 




:M''' ^ [0,M] 
for some M > 
= max(u, 0) 


-Int 


exp(— t) 


-In |^f;exp(-0,/,(x,))j 


/i,/2 iM'^^M 



of /i , /2 , . . . , /n with generating function G Q is defined as the functional 
Q^(/i,...,/„)(x) =V(/i(xi),.../„(x„)) (5) 
for x = (xi, . . . ,x„)', Xj G d = (di, . . . ,d„)' and |d| = d. 

Throughout the paper, we refer to ■i/' G H or the corresponding as 
the generating functions of Q^. Note that Qip{f, ...,/) = / for any function 
/ and generating function ip. 

Example 1. In TableUlwe present four basic examples of the quasi- arithmetic 
composition considered above. Some conventions are needed in order to deal 
with possibly ill-defined values. 

Ordering relations as well as a minimal element can be found among 
the set of quasi-arithmetic compositions of n € N fixed functions indexed 
by convex generating functions (a maximal element can also be found only 
when both fixed functions are upper bounded). From now on, we shall write 

n 

QG(.fl,...Jn)=llf!\ 

i=l 
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which is the quasi-arithmetic composition associated with ip{t) = exp(— i) 
that generates a geometric average (constraints on the weights are specified 
in Tabled]). Also, let 



for /i, ...,/„: M'^ —> [0, M]), be the quasi-arithmetic composition associated 
with ip{t) = M{1 — t/M)^ that generates the arithmetic average. Finally, 
we shall call 



the quasi-arithmetic composition associated with ip{t) = 1/t and generating 
the harmonic average. 

We shall write gi < g2 whenever 51 (xi, . . . , x„) < 52 (xi, . . . , x„) for all 
X = (xi,...,x„)' G M'^. Finally, recall that a function g is subadditive 
whenever g{a + b) < g{a) + g{b) for all a, h in its domain. Thus, we have the 
following result: 

Proposition 1. Given any set of functions fi, ■ ■ ■ , fn CLnd generating func- 
tions (/72 € $cm (which V',V'i)'02 € £}cm chfc associated with, respec- 
tively), and given the same set of weights {Oi}^^^ for any pairwise compari- 
son, we have the following point-wise order relations, 



• Qg(/i, . . . , /n) < Q^(/l, ...,/„)< Er=l (^^f^; 

• Q^fi, ■■■,fn)< Qa(/i, . . . , /n) whenever /i, ...,/„: ^ [0, M]; 

• Qvi (/ii • • • ) /n) < 2i/'2 • • • ) fn) if and only if the function Lp'^'^ o (^2 



is subadditive. 

This result extends to the functional case those reported by Nelsen 
(1999). Its proof can be found in the Appendix. 

2.3 Other useful notions and notation 

A real mapping 7 : X C M'^ ^ R is called an intrinsically stationary vari- 
ogram (Matheron, 1965) if it is conditionally negative definite, that is 



n 



QA(/l,...,/n) = 5]^i/, 



Qh(/i, ■ ■ ■ ,fn) 




n 



n 




1=1 i=j 
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for all finite collections of real weights Oj summing up to zero and all points 
Xj € X. The restriction to the isotropic case is analogical to that of covari- 
ance functions, that is 7(x) := 7(||x||), x G X and 7 : — > M conditionally 
negative definite. 

A completely monotone function if is a positive function defined on the 
positive real line and satisfying 

>0, t>0, 

for all n G N. Completely monotone functions are characterised in Bern- 
stein's theorem (see Feller, 1966, p. 439) as the Laplace transforms of positive 
and bounded measures. By a theorem of Schoenberg (1938), a function C 
is radially symmetric and positive definite on any d-dimensional Euclidean 
space M'^ if and only if C(x) := (/9(||x|p), x G M'^, with if completely mono- 
tonic on the positive real line. 

Bernstein functions are positive functions defined on the positive real 
line, whose first derivative is completely monotonic. Once again, an intimate 
connection with (negative) definiteness arises, as 7 is a radially symmetric 
and conditionally definite function on any d-dimensional Euclidean space 
if and only if 7(x) := S(||x|p), with B a Bernstein function. 

Sufficient conditions for positive definiteness are stated in Polya's criteria 
(Berg and Forst, 1975) in M^; and in Christakos (1984, 1992) and Gneiting 
(2001), who extend criteria of the Polya type to M"^. 

3 Theoretical Results 

In this section we present theoretical results in a general setting, e.g. work- 
ing with arbitrary partitions of d-dimensional spaces as explained previously. 
In particular, we shall obtain permissibility criteria for quasi- arithmetic av- 
erages of covariance functions on W^. This will be done for (a) a general 
case in which the respective arguments of the covariance functions used for 
the quasi-arithmetic average have no restrictions; (b) the restriction to the 
component-wise isotropic case; and (c) the further restriction to covariances 
that are isotropic and defined on the real line. Also, we shall show some 
properties of this construction. In particular, we refer to the associativity 
of quasi-arithmetic functionals and to the extension of ordering relations in 
Proposition [1] to the case of compositions of covariance functions. The proof 
of these new results can be found in the Appendix. 

Proposition 2. (a) General case Let G <I>cm o,nd Cj : M'^' — > M = 

1, . . . ,n) be continuous stationary covariance functions such that U"Cj(M'^*) C 
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D{i^^^) and d = (di, . . . ,dn)' , |d| = d. If the functions Xj i-^ ip^^ o 
Cj(xj), i = l,...,n, are intrinsically stationary variograms on W^*, 
then 

(Ci, . . . ,C„) (xi, . . . ,x„) (6) 
is a stationary covariance function on W^. 

(b) Component-wise isotropy Let (p G <I>cm o-n-d Cj : R"^' —>■ (i = 

1, . . . ,n) be continuous stationary and isotropic covariance functions 
such that \J'^Ci{W^^) <Z D{Lp-~^) and d = {di, . . . ,dn)' , \d\ = d. If the 
functions x i— s- 99"^ o Ci(x) are Bernstein functions on the positive real 
line, then 

Q^(Ci,...,C„)(||xi||,...,||x„||) (7) 
is a stationary and fully symmetric covariance function on . 

(c) Univariate covariances Let ip G <I>cm o-nd Cj : M — > M (i = 1, . . . ,n) 

he continuous stationary and isotropic covariance functions defined 
on the real line such that U"Ci(IR) C D[ip^^) and |d[ = n. If the 
functions x 1— s- tp^^ o Cj(x) are continuous, increasing and concave on 
the positive real line, then 

Q^(Ci,...,C„)(|xi|,...,|x„|) (8) 

is a stationary and fully symmetric covariance function on M". 

It is worth noticing that case ([8]) represents a covariance permissibihty 
condition that does not depend on the Euchdean metric, as it is function of 
the Manhattan or city block distance. For a detailed discussion of the limi- 
tations of the Euclidean norm-dependent covariances, see Banerjee (2004). 

The previous result is of particular importance, since it has implications 
both in terms of averages of covariance functions and mixtures thereof. In 
the following we derive sufficient permissibility conditions for geometric and 
harmonic averages of covariance functions. 

Corollary 1. f'Geometric average^ Let Ci : ^ "M. (i = 1, . . . ,n) be 

continuous permissible covariance functions. Let 6i (i = 1 . . . ,n) be non- 
negative weights summing up to one. If the functions x 1-^ — ln(Ci(x)), 
X > satisfy any of the relevant conditions described in (a),(b) and (c) of 
Proposition above, then 

n 

QG(ci,...,a)=ncf 

1=1 

is a covariance function. 
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An example of this setting can be found by using the function x i— > 
(1 + X positive argument, 5 G (0, 2] and e positive, also known as gen- 

eralised Cauchy class (Gneiting and Schlather, 2004). One can verify that 
the composition of this function with the natural logarithm is continuous, 
increasing and concave on the positive real line for 5 G (0, 1]. Another func- 
tion satisfying these requirements is the function x i— > exp(— x"^), which is 
completely monotonic for 5 G (0, 1]. It should be stressed that these permis- 
sibility criteria do not apply to compactly supported covariance functions, 
such as the spherical model (Christakos, 1992). 

Corollary 2. (^larmonic average^ Let Ci : R'^' ^ W (i = 1, . . . ,n) be 

continuous permissible covariance functions. If the functions x i— > Ci{x)~^ , 
X > satisfy any of the relevant conditions described in (a),(b) and (c) of 
Proposition^^ above, then for Oi >0 such that 9i = 1, the 

Qh(Ci, . . . ,C„) = 

is a covariance function. 

In order to complete the picture about quasi-arithmetic covariance func- 
tions, it would be desirable to establish at least some of their algebraic prop- 
erties. The following results show some important features of the theoretical 
construction obtained from Proposition [2l 

Proposition 3. (^Associativity^. Consider the same arrangements as in 
Proposition^^ and set 9i = 1/n (^z = 1, . . . ,n). Let ipi G <l*cm CLnd 992 £ ^• 
// the functions x 1— > ip^^ o ip2{x), x 1— > ip2^ o Cj(x) ( i = 1,... ,k) and 
X I— > ip^^ o Cj{x) (j = k + 1, . . . ,n; t > and k < n) satisfy any of the 
relevant conditions described in (a),(b) and (c) of Proposition above, 
then 

), (9) 

Qipi (Qv2(Ci, • • • , Cfc), Cfc+i, . . . , C„) (||xi||, . . . , ||xfc||, ||xfc+i||, . . . , ||x„||) 

(10) 

and 

Qvi (Qv2(Ci, • • • , Cfc), Cfc+i, . . . , C„) ...,\xk\, \xk+i\, ...A^n\) (11) 

and accordingly any coherent permutation of ipi,(p2 with Cj (i = 1, . . . ,n) 
are covariance functions. 
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One may notice that we deviate from the associativity condition defined 
in Bemporad (1926) and (in a more general form cahed decomposabihty) in 
Marichal (2000). Nevertheless, it can be readily verified that the associativ- 
ity condition of Bemporad is always satisfied for our construction, whereas 
the strong and weak decomposabihty of Marichal (2000) is certainly satisfied 
by construction but not, in general, by constructions ([6]) and ([7|). 

Proposition 4. (^Ordering and preference relations^ For any set of co- 
variance functions Ci, . . . , C„ and any arbitrary generating functions ip, (pi,(p2 G 
<I>cm with associated ilj,ipi,ip2 £ O^cm, if (i) the constraints in either (a), (b) 
or (c) of Proposition\^ are satisfied and (ii) the same set of weights 
is used for any pairwise comparison, then we have the point-wise order re- 
lations, 

(1) Qg(Ci, . . . , c„) < Q^(Ci, . . . , c„) < Er=i 

(2) Qh(Ci,...,C„) < Qa(Ci,...,C„) whenever Ci,...,Cn:M.''^ ^ [0,M]; 

(3) Q^-^(Ci, . . . ,C„) < Q^2(Ci, . . . ,C„) if and only if the function ip^^ o(f 2 
is subadditive. 

The proof is omitted as it follows the same arguments as that of Propo- 
sition [H It should be stressed that the first inequality shows that quasi- 
arithmetic covariances satisfy Jensen's inequality. This fact has some im- 
plications on the variance of the stationary random field generated by a 
quasi-arithmetic operator. This will be discussed in the subsequent sections. 
Also, it is important to specify that other results, in particular inequalities 
involving quasi-arithmetic operators, can be readily extended to the case of 
quasi-arithmetic compositions of covariance functions. This is the case, for 
instance, of Theorem 1 in Abramovich et al. (2004). 

4 Using quasi-arithmetic functionals in the con- 
struction of nonseparable space-time covariance 
functions 

4.1 Review of space-time covariance functions 

Natural (physical, health, cultural etc.) systems involve various attributes, 
such as atmospheric pollutant concentrations, precipitation fields, income 
distributions, and mortality fields. These attributes are characterised by 
spatial-temporal variability and uncertainty that may be due to epistemic 
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and ontologic factors. In view of the prohibiting costs of spatially dense 
monitoring networks, one often aims to develop a mathematical model of 
the natural system in a continuous space-time domain, based on sequen- 
tial observations at a limited number of monitoring stations. This kind of 
problem has been a motivation for the development of the spatio-temporal 
random field {S/TRF) theory; see Christakos (1990, 1991, 1992) and Chris- 
takos and Hristopulos (1998) for a detailed discussion of the ordinary and 
generalised S/TRF theory and its various applications. In the following, 
we slightly change notation in order to be consistent with classical nomen- 
clature in the Geostatistical literature (Cressie and Huang, 1999; Gneiting, 
2002). Let {Z(s,t), (s,i) G M'^ x M} be a real-valued S/TRF, where s,t 
denote respectively the spatial and temporal position. Then, the function 
Cs^t(si, ti, S2, ^2) = cov(Z(si, ti), Z(s2, ^2)), defined on the product space 
M'^ X M X M'^ X R, is the covariance function of the associated S/TRF Z if and 
only if it is permissible, i.e. satisfies ([T]). When referring to the spatial index, 
the term homogeneity instead of weak stationarity is equivalently adopted 
(Christakos, 1990, 1992). Thus, under the assumption of spatial homogene- 
ity and temporal stationarity (sometimes, simply called spatio-temporal sta- 
tionarity in the weak sense), the underlying S/TRF, ergo denoted by H/S, 
has finite and constant mean and the covariance function, defined on the 
product space R'^xR, is such that cov(Z(si, fi), Z(s2, t2)) = Cs^j(h, u), with 
(h, u) = (sj — Sj,ti — tj) E M'^ X R denoting the spatio-temporal separation 
vector, and the h and u denoting the spatial and the temporal lags, respec- 
tively. The special case of isotropy in the spatial component and symmetry 
in the the temporal one is denoted as 

C,,t(h,u) :=a,t(||h||,|n|), (12) 

where the ||.|| denotes the Euclidean norm. Obviously, (jl2p is fully symmet- 
ric. 

Another popular assumptions concerning the S/TRF model is that of 
separability, that is (Christakos and Hristopulos, 1998; Mitchell et al, 2004) 

C,,t(h,u) = C,,t(h,0)C,,t(O,^i). (13) 

In other words, separability means that the spatio-temporal covariance struc- 
ture factors into a purely spatial and a purely temporal component, which 
allows for computationally efficient estimation and inference. Consequently, 
separable covariance models have been popular even in situations in which 
they are not physically justifiable. Another interesting aspect is that sep- 
arable covariances are also fully symmetric, whereas the converse is not 
necessarily true. 
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It has been argued in the relevant hterature that separable models allow 
for ease of computation and dimensionality reduction, as the space-time co- 
variance matrix is obtained through the Kronecker product of the marginal 
spatial and temporal ones. However, separability is an unrealistic assump- 
tion for many applications, since it implies a considerable loss of information 
about important interactions between the spatial and temporal variations. 
Therefore, various techniques have been introduced for generating different 
classes of nonseparable spatio-temporal covariance models. Most of these 
techniques have been developed in the context of applied stochastics anal- 
ysis and include H/S as well as non-H/S covariances (e.g., Yaglom, 1948; 
Gandin and Boltenkov, 1967; Monin and Yaglom, 1965, 67; Christakos, 
1990, 1991, 1992, 2000; Christakos and Hristopulos, 1998; Christakos et al, 
2002, 2005; Kolovos et al, 2004). Also, covariance models have been de- 
veloped in the context of spatio-temporal statistics (e.g., Jones and Zhang, 
1997; Cressie and Huang, 1999; De Cesare et al, 2001; Fucntes and Smith, 
2001; Gneiting, 2002; Ma, 2003; Fernandez-Casal et al, 2003; Stein, 2005a; 
Porcu et al, 2006a). 

Being the Laplace transform of positive bounded measures, completely 
monotone functions are particularly appealing for the construction of sta- 
tionary and nonstationary space-time covariances. In particular, they are 
intimately connected with the concept of mixture-based covariance functions, 
that has been repeatedly used by several authors. In the stationary case, 
see Gneiting (2002), Ma (2002;2003), Ehm et al. (2003), Fernandez-Casal 
et al. (2003), Porcu et al. (2006a). In the nonstationary case, Paciorek and 
Schervish (2004;2006), Stein (2005b), and Porcu et al. (2006b) have made 
use of this technique. 

Also, mixture-based covariances have been developed with less sophis- 
ticated instrTinicnts than completely monotone functions. This is the case 
of the so-called product sum model (De Cesare, 2001) and their extensions 
(De laco et al, 2002; Porcu et al, 2006c). This group of authors builds 
nonseparable space-time covariances through simple application of the basic 
properties of covariances seen as a convex cone. The mixture-based proce- 
dures and the basic properties of covariance functions are properly combined 
in the present paper. 

In the following, a stationary RF with a quasi-arithmetic covariance 
function will be called quasi- arithmetic random field and denoted with the 
acronym QARF. 
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4.2 On the representation and smoothness properties of QARF 

In this section we focus on the representation of a QARF and discuss the 
smoothness properties in terms of (mean square) partial differentiability. 
These properties are intimately related to those of the associated covariance 
function. 

Let Zi{s) be univariate mutually independent continuous weakly station- 
ary Gaussian random processes defined on the real line {i = 1, . . . , (d + 1); 

s G R; and d € In particular, let the process Z^^i be continuously 

indexed by time t. Consider also a {d + l)-dimensional nonnegative random 
vector R = {Ri, . . . , Rd+i) with Ri independent of Zi. Let the univariate 
covariances C^. and the temporal covariance Ct be respectively associated 
with Zi, i = 1, . . . ,d and Zd+i- In the following we shall assume these covari- 
ances to be stationary, symmetric, and of the type Csi{hi) = exp(— 
i = 1, . . . ,d, and Ct{u) = exp(— with Vi = ip~^ o Q, Vf = ° Q, 
where the (p G $cm and Cj are positive definite and such that the compo- 
sitions Vi are continuoTis, increasing and concave on the positive real line. 
Positive definiteness of this construction is guaranteed by direct application 
of the theorem of Schoenberg (1938) and according to a Polya type criterion 
(see Berg and Forst (1975), Proposition 10.6). 

We are interested in inspecting the properties of the following stationary 
spatio-temporal scale mixture-based random field, defined on M*^ x R, 



Z{s,t) = Zd+i{tRd+i)f[Zi{siRi), (14) 



1=1 



with s = {si, . . . , Sd)' G R"* and i G R. It can be easily seen that the 
covariance structure associated to this random field is nonseparable, as 

Cs,t{h,u) = [ exp ( -Y^UidhiDn - ut{\u\)rd+?\ dF{r), (15) 

with h = {hi,...,hd)' G R'^, u G R and r = (n, . . . , r^+i)' G R''+\ and 
where F is the {d + l)-variate distribution function associated to the ran- 
dom vector R. If F is absolutely continuous with respect to the Lebesgue 
measure, then previous representation can be reformulated with respect to 
the {d + l)-variate density, say /, that is 

Cs,t(h,u) = / cxp - Vi/id/iiDn - i^t{\u\)rd+i I f{r)dr. 

M^' V i=i J 
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It can be seen easily that this construction ahows for the case of sep- 
arabihty if and only if the integrating (d + l)-dimensional measure F (or 
equivalently its associated density /) factorises into the product of {d + 1) 
marginal ones, i.e. if the nonnegative random vector R has mutually inde- 
pendent components. 

Now, if we suppose that the measure F is concentrated on the line ri = 
... = Vd+i = r and that ^p'"^ is such that := C[F], i.e. the Laplace 
transform of the positive univariate measure -F, then we obtain that the 
QARF is a special case of (jlSp . For this random field, further inferences 
may be made about its mean square partial differentiability if, additionally, 
we assume the function t h-> exp(— z/i(t)) to be absolutely integrable on the 
positive real line (i = 1, . . . , In this case, one can show (technicalities 

can be found in the appendix) that the fe-th order mean square partial 
derivative with respect to the i-th coordinate exists and is finite whenever 
the function xi^) = J^i^^^-^u;f''ci{oJi;r)du;i, with Cj := .F"^[exp(-rfj)] = 
Jjg exp(— iwj/ij — rvi{hi))Ahi is measurable with respect to F (i = 1, . . . , d). 

It should be stressed that Proposition S] gives some more information 
about the characteristics of the underlying QARF in terms of variance, as 
QARF can be ordered with respect to their minimum or maximum variance. 
Thus, the QARF generated hy ip := A has the largest variance among all 
the other QARF, for any choice of ip. 

Finally, observe that the trivial quasi-arithmetic composition, obtained 
by setting Cj := 99 G ^>cm (i = 1, . . . , d) preserves permissibility and results 
in a model of the type 

C,(h) = ^ (0'h) 

where = {9i, . . . , 9d)' G IK'^ and h = (hi, . . . , h^)' is the arbitrary partition 
of the spatial lag vector. Thus, the trivial case results in a composition of a 
completely monotonic function with an affine function, and can be used for 
modelling geometric anisotropics. 

4.3 Applications of quasi-arithmeticity in the construction 
of stationary nonseparable space-time covariances 

The results presented in the previous section can be useful in the construc- 
tion of space-time covariance functions. For this purpose, some considera- 
tions are in order. The functional in equation ([4]) should be adapted in the 
spatio-temporal case. It does not make sense to consider a weighted average 
of a spatial covariance with a temporal one. Also, the use of weights forbids, 
in this case, the construction of nonseparable model admitting separability 
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as a special case. For this reason, we suggest to reduce the class in equation 
dH) to the case of trivial weights, e.g. 9i = 1, Vi. Then, it is more appropriate 
to call the class Archimedean, in analogy with the class built in Genest 
and MacKay (1986). It should be mentioned that one can easily prove that 
the restriction to trivial weights does not affect the permissibility of the 
resulting covariance function, provided that either one of the constraints 
imposed in cases (a), (b) or (c) of Proposition [2] are fulfilled. By analogy 
with the construction of copulas, and following Genest and MacKay (1986), 
the definition of generator for (p G $cm is more meaningful. 

Working this way, one can obtain some new families of covariance func- 
tions whose analytical expressions are familiar for those interested in proba- 
bilistic modelling through copulas. It should be noticed that all the families 
we propose in this section include separable covariance models as special 
cases, depending on the parameter values of the generators. 

Example 1. The Clayton family 

Consider the completely monotone function 

<f{x) = (l + xy^/^', x>0, (16) 

where Ai is a nonnegative parameter with inverse ip^^{y) = y~^^ — 1. Ob- 
serve that ()16p is the generator of the Clayton family of copulas; refer to 
Genest and MacKay (1986) for mathematical details about this class. 

Also, consider the covariance functions Cs(h) = (1 + ||h||)-VA2 and 
Ct(«) = (1 + |n|)~"^/'^3, for A2,A3 positive parameters. It is easy to ver- 
ify that ip~^ (Cj(y)) = (1 + y)'^^l^i {i = 2,3) is a positive function whose 
first derivative is completely monotone if and only if Ai < Aj. Under this 
constraint, and applying case (b) of Proposition [21 we find that 

C,,t(h, u) = Q^(C„ Ct)(h, u) = a" [(1 + \HT + (1 + W " 1]"'^^' 

(17) 

is a valid nonseparable stationary fully symmetric spatio-temporal covari- 
ance function, with pi = Ai/Aj, Aj > (i = 2,3) and o"^ a nonnegative 
parameter denoting the variance of the underlying S/T process. It is inter- 
esting that both margins (the spatial and the temporal one) belong to the 
generalised Cauchy class. This is desirable for those interested in the local 
and global behaviour of the underlying random field. 

Another covariance function that preserves Cauchy margins can be ob- 
tained through the following iter. Consider the function x i— s- ip{x) = x~°', 
t > 0, that belongs to <&cm for any positive a, being the Laplace transform 
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of the function x ^ a;"~^/r(a), with a positive parameter and V the Euler 
Gamma function. Also, for the spatial margin consider Cs(h) = (1+ ||h||^)~'^ 
that belongs to if and only if 5 € (0, 2] and e is strictly positive. Finally, 
let the temporal margin be of the type Ct{u) = \u\^f, which is not a station- 
ary covariance function but respects the composition criteria, as ip~^ o Ct 
is continuous, increasing and concave on the positive real line if and only if 
a < p. Then, it is easy to prove that Q^{Cs, Ct)(h, u) is a valid space-time 
function if and only if a < e and a < p, and that this covariance function 
has margins of the Cauchy type. 

Example 2. The Gumbel-Hougard family 

Consider the completely monotone function, the so-called positive stable 
Laplace transform 

^(x) =e-^'^^S x>0, (18) 

where Ai > 1. Equation (18) admits the inverse ip'^iy) = {—\n.{y))^'^ . This 
is the generator of the Gumbel-Hougard family of copulas, whose mathemat- 
ical construction and details are exhaustively described in Nelsen (1999). 
Consider two respectively spatial and temporal covariance functions admit- 
ting the same analytical form, i.e. Cs(h) = exp(— ||h||-^/'^2) ^q^j Ct(n) = 
exp(— l^l^/'^''*). Then, it can be easily verified that ip~^ {Cs{y)) = y^^/^^, 
and Lp~^ (Ct(y)) = y^^^-*"'-^ are always positive for y > and possess com- 
pletely monotone first derivatives if and only if Ai < Aj, i = 2, 3. So we get 
that 

Cs,t{h,u) = Q4Cs,Ct){h,u) = a^exp (-(Hhfi + \ur)'/^') , (19) 

is a permissible nonseparable stationary fully symmetric spatio-temporal 
covariance function, with pi = Xi/Xi, i = 2,3 and as before. 

Example 3. The power series family 

The so-called power series 

ip{x) = 1 - (1 - exp{-x))^/^\ x>0, (20) 

with Ai > 1, admits the inverse ip~^{y) = — ln(l — (1 — y)'^^)- Suppose a 
spatial and a temporal covariance function of the same analytical form as 
in Example 2. The composition </?~^(Ci(y)) = — ln(l — (1 — exp(— y))^^/^*). 
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i = 2,3, is always positive for y > and admits a completely monotone first 
derivative if and only if Ai < Aj. So we get that 

Qv>(Cs,CT)(h,u) = l-(l-exp(-||h||)r -(l-exp(-l^xl)r 

+(1 - exp(-||h||))^i(l - exp(-|n|)r (21) 

is a permissible nonseparable stationary fully symmetric spatio-temporal 
covariance function, with pi = Ai/Aj (i = 2, 3). 

Example 4. The semiparametric Frank family 

Here we show that a nonstationary covariance function can be obtained, 
starting from the proposed setting, even if the arguments of the quasi- 
arithmetic functional are not covariance functions. 

The Frank family of copulas (Genest, 1987) is generated by the function 
(p{x) = jln (l — (1 — e^'^)e^^) with inverse (f^^iy) = —In ((1 — e~'^^) /(I — e^ 
Nelsen (1999) shows that, for A positive, (p is the composition of an abso- 
lutely monotonic function with a completely monotonic one, i.e. a com- 
pletely monotonic function. As far as the inverse is concerned, it is easy to 
see that (^^^07 (7 an intrinsically stationary variogram) is negative definite. 
Thus, the 

(1 _e-^7s(h))(i _e-^^TW)\ 



for 75, 7t intrinsically stationary variograms defined on and M, respec- 
tively, is a stationary space-time covariance function. 

4.4 Quasi-arithmeticity and nonstationarity in space 

The direct construction of spatial covarianccs that arc nonstationary is any- 
thing but a trivial fact. Only few contributions refer to this kind of construc- 
tion: see, among them, Christakos (1990, 1991), Christakos and Hristopou- 
los (1998) and Kolovos et al. (2004). More recent contributions can be 
found in Stein (2005b) and Paciorck and Schcrvish (2006). It seems that 
something more could be done concerning the construction of nonstationary 
covariances, knowing that stationarity is very often an unrealistic assump- 
tion for several physical and natural processes. 

In this section, we show that both approaches discussed in Stein (2005b) 
and Paciorek and Schervish (2006), admit a natural extension through the 



Cs,t(h, u) 



-In 



21 



use of quasi-arithmetic functionals. We need to introduce some more no- 
tation concerning the restriction of the class ^I'cm- For abuse of notation, 
a completely monotonic function is the Laplace transform of some positive 
and bounded measure, so that ip G $cm if and only \i ip := C[F]. 

Theorem 1. Let he a mapping from W x W to positive definite p x p 
matrices, F a nonnegative measure on M+, S $cm and g a nonnegative 

function such that, for any fixed s G W, = 923^ ° dU^) ^ ^H^)- Define 
E(si,S2) = l/2(S(si) + S(s2)) andQ{si,S2) = (si - S2)'S(si, sa)^^^! - S2). 
Then, 

a(si,S2) = T'/'i (Q(si,S2)r) Q^,(5s„fe)(r)dF(r) 

|S(Si,S2)|V^ Jq 

(22) 

with 0i = 1, i = 1,2, is a nonstationary covariance function on M.'^ x M^. 

Some comments are in order. One can see that Stein's (2005b) result 
is a special case of (p2]) . under the choice ipi{t) = exp(— t), t positive, and 
ip2 ■= G. So is the result in Paciorek and Schervish (2006). Nevertheless, the 
form we propose has some drawbacks that need to be noticed explicitly. The 
first problem is that it is very difficult to obtain a closed form for expression 
(j22|) . unless one chooses ip2 ■= G. Secondly, if one's purpose is to generalise 
the Matern covariance function (Matern, 1960) to the nonstationary case, 
then the problem has already been solved by Stein (2005b) and his approach 
is in our opinion the most suitable, as he finds a Matern type covariance that 
allows for a spatially varying smoothness parameter and for local geomet- 
ric anisotropy. On the one hand, the Matern covariance possesses some 
desirable features (i.e., it allows for arbitrary levels of smoothness of the 
associated random field). On the other hand, there are other covariance 
functions that are of considerable interest to the statistical and scientific 
communities, such as the nonstationary ones. In particular, we refer to the 
Cauchy class, whose properties (in terms of decoupling of the global and local 
behaviour of the associated random field) have been discussed in Gneiting 
and Schlather (2004). After several trials, we did not succeed in obtaining a 
Cauchy type nonstationary covariance. Several examples of covariances can 
be derived that belong to the class (j22l) and can be numerically integrated. 
Here, we propose some closed form that is obtained by letting ^2 '■= G. 
As a first example, take d-F(r) = dr, 921 (t) = t^~^, that is completely 
monotonic for A G (0, 1), g{T; ai, Vi) = (1 + a(sj)T)^'^(^»\ where a and u are 
supposed to be strictly positive functions of Sj (i = 1,2) and additionally 
< a(sj), i^(sj) < TT. One can readily verify that all integrability conditions 
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in Theorem [T] are satisfied. Letting k = ^^^^^g'^ s^)ji/2 ^ Q(si,S2)^ ^ and 
using [3.259.3] in Gradshteyn and Ryzhik (1980), one obtains the following 
covariance function 

C,(si, S2) = A:a(si)-^B(A, u{si)+u{s2)-X)2Fi (u{s2), A; u{si) + 1/(83); 1 - 4^ 

where B(.,.) is the Beta function, and 2Fi{-, ■,-,■) is the Gauss hypergeo- 
metric function. 

Another example is obtained by considering F(dr) = exp(— r)dr, (Pi{t) = 
t"~^i with u G (0,1), g{T;Si) = exp(— ^^^r), with a(sj) strictly positive 
(i = 1, 2) and using [3.478.4] of Gradshteyn and Ryzhik (1980), we find that 

c.(s„ s,) = 2k (5!M±i^) "^'^ x:. (^2 (5!(!i)±i^) '^'^ 

is a nonstationary spatial covariance that allows for local geometric anisotropy, 
but has a fixed smoothing parameter. 

It should be stressed that numerical integration under the setting (j22]) 
could outperform previously proposed models if the objective is to find a dif- 
ferent type of interaction between the local parameters characterising the in- 
tegrating function g. In all the examples proposed by Paciorek and Schervish 
(2006) and Stein (2005b) the varying-smoothing, spatially adaptive param- 
eter is obtained as the semi-sum of a parameter acting on the location si 
with another one depending on the location 82- This is a serious limitation 
of the method, as only one type of interaction can be achieved. Starting 
from a different setting, Pintore and Holmes (2004) obtain the same type 
of spatially adaptive smoothing parameters. Thus, quasi-arithmetic func- 
tionals could be of help, at least through numerical integration. Finding a 
closed form for a functional different than Qq remains an open problem to 
which we did not find any solution for the meantime. 

5 Conclusions and discussion 

In this work, novel results are presented concerning permissible spatial- 
temporal covariance functions in terms of the theory of quasi arithmetic 
means, and valuable insight is gained about their space-time structure. The 
theory of quasi arithmetic means is used, together with the relevant permis- 
sibility criteria, to derive new classes of nonseparable space-time covariances 
and to investigate their properties in considerable detail. 
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There are several possible avenues for research based on the results of 
the present paper. From the spatial and spatiotemporal statistics perspec- 
tives, the QARF representation considered in this paper seems very promis- 
ing. E.g., the QARF may provide a naive separability testing procedure, 
as follows. Consider the set of all possible generators of the QARF 
class and let Q be the set of parameter vectors indexing the generators, i.e. 
Q : {6 £ W : ifQ G ^*cm}- Thus, testing for separability of the covariance 
function associated to the QARF is equivalent to testing for the null hypoth- 
esis ifo := ipo{x]6Q) = exp{—dQx), x > 0. This is a topic worth of further 
investigation. 

Another interesting topic is that of generator estimation. A tempting 
choice would be to consider the approach proposed by Genest and Rivest 
(1993), and that in Ferguson et al. (2000) who find a serial version of the 
Kendall's tau and a relationship between this concordance index and the 
generating function ip. In this case, one would extend the result of these 
authors at least to the lattice 7?. 

Future research effort would also focus on a deeper study of the prop- 
erties of quasi-arithmetic covariances along several directions. For instance, 
it would be interesting to find limit properties of the covariance functions 
specified through this construction. 
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6 Appendix 

6.1 Proof of Proposition [1] 

Let us begin by showing the final part, i.e. that 



Qipi ifl, - ■ ■ ,fn) < (/l) • • • 1 fn) 
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for V'l ) '02 £ Hem ) if and only if if^ 
Let g = fi^ o f2 and denote Xj 
Thus, (123]) is equivalent to 



ip2 is subadditive 

ip^-^ O /i(Xi). 




(24) 
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For the necessity, assume Aip^ < A^p^- Applying tp^ in both sides of (p^ 
we obtain the result. 

For the if part, assume holds. Thus, applying ipi to both sides of 
the same inequality the result holds. 

Let us now prove that ^n(/i, •••,/«)< A^{fi, ...,/„)< 9ifi. 

The lower bound in the inequality is direct consequence of Lemma 4.4.3, 
Corollary 4.6.3 in Nelsen (1999, p. 110) and references therein, while the 
upper bound is a direct property of convex functions. 

6.2 Proof of Proposition 2 

Recall that, for Bernstein's theorem (Feller, 1966, p. 439), f € $cm if and 
only if 



Jo 

with F a positive and bounded measure. Now, for the proof of (a), notice 
that, if 99 G ^cm, equation Q can be written as 



Now, observe that if, for every i, ip o Cj is a variogram, then so is the 
sum Y17=i^i^~^ ° ^i- Then, by Schoenberg theorem (1939), we have that 
for every positive r, the integrand in the formula above is a covariance 
function. So is the positive scale mixture of covariances. This completes the 
proof. 

For (b), it is sufficient to notice that, by Schoenberg theorem (1938), 
the mapping 7 : R"' ^ is the variogram associated to an intrinsically 
stationary and isotropic random field if and only if 7(h) = ■i/'(||h|p), for a 
Bernstein function. The rest of the proof comes from the same arguments 
of point (a). 

For (c), notice that, being continuous, increasing and concave on [0, 00), 
each of the functions ip~^{Ci), is negative definite on M, according to a Polya 
type criterion (see Berg and Forst (1975), Proposition 10.6). So is their sum. 
Since (p G ^cnu by Schoenberg's theorem (cf. Berg and Forst, 1975), we get 
once again the result. 

6.3 Proof of Proposition [3] 

We shall only prove the result in ([9]), as (jlOp and (jlip follow the same argu- 
ment. Recall that here we impose 6i = 1/nVi. Canh(i) = (||hi ||, . . . , ||hfc||)' 




(25) 
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and h(2) = (||hfc+i||, . . . , ||h„||)', k<n. Also, let /i = l/nY:]=k+i Vi^ ° Cj, 

/2 = ° V2 and /a = 1/?^ X]i=i ^ ° Clj. Thus, equation ([9]) can be 
written, using Bernstein's theorem, as 

poo 

/ e-r/l(h(2))-r/2o/3(h(l))^_p_^(^^) 

where -Fi is the distribution associated to its Laplace transform ipi. Thus, 
the proof follows straight by arguments of the previous proofs. 



6.4 Proof of Theorem [T] 

The result is a direct consequence of the work of Stein (2005b). First, 
observe that C(si, si) = Jg g{T; si)dF{si) < 00. Now, we need to show that 
"^^=1 aiajC{si,Sj) > 0, for every finite system of arbitrary real constants a,, 
i = 1, . . . ,n. Write Kj'^'^ for the normal density centered at Sj and with 
covariance matrix (rrr')~^Ssj. Also, by abuse of notation, := C[Fi] and 
ip2 := >C[F2]. Then, by a convolution argument in Paciorek (2003, p. 27), by 
Bernstein's theorem, and using Fubini's theorem, we get 

i,s = l 

= y a,a, IE,. I i IE3. 1 i P P P ^-■•■•'■'(u)if-''--'-'{u)du') (rrr')-''/2o-'-'<l''"'^2~'°s(-'«i) + l/2»'2"'°9(-i=j»dF(r)di 

j"^! ' ^ Jo Jo Jo \JKP ' ' J 

= n"^^ [°° [°° [°° [ (ai|Es.|iK^'''''^'{u)o"'''^/^^2"^°»<^'°i>Vdu(Trr')"''/^dF(T)dFi{r)dF2(r') > 
Jo Jo Jo JRP \ ' ' / 

Thus, the proof is completed. 



6.5 M.S. differentiability of the QARF 

Recall that the existence of the k-th order i-th mean square partial derivative 
of Z is related to the existence of the 2fc-th order i-th mean square partial 
derivative of the covariance function. This is in turn related to the spectral 
moments through the following formula (Adler, 1981, p. 31): 

(-l)'^^lh=o = / -I'dCiu.) < 00. 

Under the setting in Section 4.2, we ensure the existence of the spectral 
density c, thus 
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/ a;f c(a;)da; oc / a;f / e-*'^''^C(h)dhdu; 

/•oo 

^'h / e-^^i^''^*(l'^*l)dF(r)dhd4^ 
Jo 



2k I „-ia)'l 



oo 

, 2k 



where h = (/ii, . . . , ^j+i, . . . , /i^)' G M*^"^ and where we repeatedly 
make use of Fubini's theorem for standard integrabihty criteria. Now, ob- 
serve that, for the assumption of absolute integrabihty of exp(— r^ifi(a;)) for 
any positive r, the last equality can be written as 




a;^- c(u;; r)c{uji; r)du;dF{r) 



where c(cj;r) := JT ^[Cr](a>), with Cr{h) = exp{—rJ2j_^-6ji'j[\hj\)) and 
where c{oJi;r) = J^~^[Cr]{oJi), with Cr{hi) = exp(— r^ii/j(|/ij|)). Then, by 
noticing that Jj^a-i c(u»;r)du; = 0^(0) < oo, one obtains that the integral 
above is finite if and only if the function 

r I— > Jj^u;^^c{ujf,r)duji is F-measurable. Thus, the proof is completed. 
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